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Abstract 


A numerical method for solving Fredholm and Volterra integral equations 
of the second kind is presented. The method is based on the use of the 
Newton—Cotes quadrature rule and Lagrange interpolation polynomials. 
By the proposed method, the main problem is reduced to solve some non- 
linear algebraic equations that can be solved by Newton’s method. Also, 
we prove some statements about the convergence of the method. It is 
shown that the approximated solution is uniformly convergent to the exact 
solution. In addition, to demonstrate the efficiency and applicability of the 
proposed method, several numerical examples are included, which confirms 
the convergence results. 
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1 Introduction 


Integral equations have lots of applications in science and engineering. Fred- 
holm integral equations arising in the theory of signal processing, which is a 
subfield of mathematics, information, and electrical engineering. In physics, 
the solution of such integral equations allows for experimental spectra to be 
related to various underlying distributions, for example, the mass distribu- 
tion of polymers in a polymeric melt. They also prevalently appear in linear 
forward modeling and inverse problems. Also, Volterra integral equations 
arise in many scientific applications such as the population dynamics, the 
spread of epidemics, and semi-conductor devices. It was also shown that 
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Volterra integral equations can be derived from initial value problems; see 
[21]. Therefore, providing effective methods for solving such equations is 
consequential. Fredholm and Volterra integral equations have been exten- 
sively studied in many papers, and numerical methods have been widely 
used to solve such equations, for instance, Newton—Kantorovich and Haar 
wavelets for nonlinear Fredholm and Volterra integral equations [1], hat ba- 
sis functions for the system of linear and nonlinear integral equations [2], 
Haar wavelets for nonlinear Fredholm integral equations [3], alternative Leg- 
endre polynomials for nonlinear Volterra-Hammerstein integral equations [5], 
operational matrix approach for nonlinear Volterra-Fredholm integral equa- 
tions [4], cubic spline method for Fredholm integral equations [6], implicitly 
linear collocation method for nonlinear Volterra equations [7], Chebyshev ap- 
proximation for nonlinear Fredholm-Volterra integral equations [9], Cheby- 
shev collocation method for the class of Fredholm integral equations with 
highly oscillatory kernels [10], collocation-type method for Hammerstein in- 
tegral equations [12], Sinc-collocation method for nonlinear Fredholm integral 
equations with weakly singular kernel [15], triangular functions for nonlin- 
ear Volterra-Fredholm integral equations [13], wavelets-Galerkin method and 
wavelets precondition for first kind Fredholm integral equations [14], spec- 
tral collocation method for Fredholm integral equations on the half-line [16], 
hybrid Legendre Block-Pulse functions for the system of nonlinear Fredholm— 
Hammerstein integral equations [17], single-term Walsh series for nonlinear 
Volterra-Hammerstein integral equations [18], piecewise constant functions 
method for nonlinear Fredholm-Volterra integral equations [19], and so on. 
The existence and uniqueness of such equations were discussed in [8, 11]. 

In this work, we consider the Fredholm and Volterra integral equations of 
Hammerstein type 


b 
u(x) = v(a) +f K(x, t)w(t, u(t))dt, 2 € {a, BI, 


u(a) = v(x) + i. K(x, t)w(t, u(t))dt, 2x € [a,b], 


where v € L?/a,b] and « € L?{a,b]? are known functions, ~ is a given non- 
linear function defined on [a, 6], and u is unknown to be determined. 


2 Method of solution 


In this section, first, we describe the Newton—Cotes quadrature rule. 
Let the interval [a,b] be partitioned into n sub-intervals, by n equally 
spaced points. That is, 


t=a, =X +th, fori=0,1,...,n, (1) 
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In—-X _ b-a 
n TS 


where the step size h is defined by h = 


The Newton—Cotes quadrature rule for a function f defined on [a, b] with 


known values at equally spaced points x;, i = 0,1,...,n, is as follows: 
b n 
[ f@dex Dwiftes), (2) 
a i=0 


for the set of weights {w;}?_). Now consider the Lagrange basis polynomials 


i;(t) = I 253 itt) and let p,(t) be the interpolation polynomial in the 


Lagrange form for the given data points (to, f(to)), (1, f(t1)),---, (tr; f(tn))- 


Then 
[soa f ontoa 
-[ (> ren) dt 
: 
= Dats f 1;(t)dt 
= Yast) (3) 
where w; = f” l;(t)dt. 


Note that we take the Lagrange interpolation points t; to be the same 
as the points x; for the Newton—Cotes quadrature rule. Now consider the 
following Fredholm and Volterra integral equations of the second kind 


b 
u(x) = (a) +f K(x, t)w(t, u(t))dt, x € [a,b], (4) 
u(x) = v(x) + A K(x, t)w(t, u(t))dt, a € [a,b]. (5) 
These can be written as 
b 
u(x) = v(x) +/ K(x, t)W(t)dt, (6) 
u(x) = v(x) + if K(x, t) UV (t)dt, (7) 


where W(t) = w(t, u(t)). By considering u,,(a) as an approximation for u(x), 
we can turn equations (6) and (7) into the following equations: 
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ico wai i * st) Uy (t)at (8) 
tin (a) = v(0) + | * (at) Vn (t)dt, (9) 
where V(t) = V(t, un (t)), which immediately implies 
v(t) =¥ (10 [/stwatnnr) (10) 
v,(t) = (« v(t) + [« (t,2)W,(a Jar). (11) 


Using quadrature formula (2) to evaluate the integral part of (10), we obtain 


b n 
/ K(t,2)V,,(a)dx = S wit(t, 4), (24) 
a 1=0 
= So win(t, Li) Uni (12) 
i=0 


where U,,; = U,(x;) and uw; = Hi I,(t)dt. Similarly, using the Lagrange in- 
terpolation for the integrand of (11) gives 


i ih, ade= to = [ Sxtiay) Paiute 


i=0 
=yo8 a I,(x)dax 
= S\wi(t)a(t, xj) Vpn (xi) 


wi (t)K(t, 21) Una, (13) 
1=0 


where w;(t =). l(x)dx and W,,; = W,,(x;). Substituting (12) into (10) and 
(13) into 6 we ee 


U,(t) =a ( v(t) + Dwin(t “).) (14) 
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respectively. Evaluating (14) and (15) at the points t = 2;, 7 = 0,1,...,n 
(the points for Newton—Cotes quadrature rule), respectively, gives 


Wnj =Ww (2000) + So santose ns) » j=0,1,...,2, (16) 


i=0 


and 


Ug = (- Ce Sane] , jg=9,1,...,n, (17) 


i=0 


where w;,; = w;(aj) = J” li(x)dx. Nonlinear systems of algebraic equations 
(16) and (17) can be solved by numerical methods such as Newton’s method. 
By solving the above systems, the values U,,;, 2 = 0,1,...,, will be known. 
Finally, by substituting (12) into (8) and (13) into (9), we find the numerical 
solutions of the integral equations (4) and (5) by 


n 


Un(x) = v(x) + So win(2, ti) Un (18) 


i=0 
and : 
Un (2) = v(x) + S° wi(x)K(a, ti) Vn i, (19) 
i=0 
respectively, where t; = 4;,i1 =0,1,...,n. 


3 Convergence of the method 


In this section, we analyze the convergence of the prescribed method in Sec- 
tion 2, which enables us to control the estimated errors. First, we provide an 
interpolation polynomial error bound, which is given in the following theo- 
rem. 


Theorem 1. [20] Suppose that f € C"*!{a,6], and let p, be a polyno- 
mial of degree < n that interpolates the function f at n+ 1 distinct points 
X0,21,---,2n € [a,b]. Then for each x € [a,b], there exists a point ¢, € [a,b] 
such that 

7(a) 


Fle) — Pal) = ey fh M (Ge), (20) 


where m(x) = I%9(a — x;). 


Theorem 2. In the case of equally spaced interpolation points zr) = a and 
Xu; = 2% + th, fori =0,1,...n, where h = baa we have 
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Proof. Suppose x € [a;,2i41]. Then for the first 7 terms of 7(a), that is 
ID 6(@ — £;) = (@ — xo)(@ — @1)...(@ — a1), and due to equally spaced 
points x;, we have 
|a — @-3| < viga — Ui; 
=(j+1)h, j7=1,2,...,4. (22) 
Thus 
MS (a — a5) < G+ DIA, (23) 
For the next two terms of z(x), that is (« — 2;)(a — v;41), and using the 


2 
simple identity a6 < (244) , we can write 


|(a — x5) (@ — wi41)| = (@ — 24) (2i41 — 2) (aj << 2441) 
= ( ore :) 
= a (24) 


For the n — i — 1 remaining terms of (zx), that is U}_j4.(% — 23) = (@ — 
Xj4+2)(@ — 2j43)...(@— &,), we may proceed as follows: 


|x = Li+;| = Vi4j —2 (x; <a“ 541) 


< Vj4;5 — Lo 


=(it+j)h, j =2,3,...,n—-i. (25) 
Thus 
n n! n—i-1 

[Waipo(e — 23)| < G+)! h . (26) 

Therefore combining (23), (24), and (26) leads to 

|r(x)| = [IIo (@ — x4)| 

! 

z qe (27) 


Theorem 3. Suppose that « € O"T'{a, b]? and that w in (4) is a function in 
C"*1la, b] with n > 0. If u(«), the exact solution, and u,,(a), the approximate 
solution defined by (18), are both in C"*1[a,6], then un(x) is uniformly 
convergent to u(). 


Proof. From (6) and (18), for every x € [a, b], we have 
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u(a) — Un (x) = ii 


- / , (seme > > rae) dt, (28) 


i=0 


b 


n b 
K(x, t)U(t)dt — So wiK(2, ti) Une (wi = / 1;(t)dt) 
i=0 a 


but So") K(2,ti) Wn ili(t) interpolates K(x,t)U(t) at the points tj, 7 = 
0,1,...,n, (¢; = 2;, the points for Newton—Cotes quadrature rule). If we 
set F(a, t) = «(a,t)U(t), from (28) and Theorem 1, then 


b x n+1 
u(a) ~un(z) = [ (z Oar e.G)) dt, ¢,€ [a,b], (29) 


where m(t) = II?_9(t — t;). Since ~,u € C”*![a,b], there is some M, > 
0 with |W("+)(t)| < M, for all t € [a,b], where U(t) = Y(t, u(t)). Also 


k € C"*1[a, b]?; thus there is some Mz > 0 with | are (a, t)| < Mpg for all 


x,t € [a,b]. Then F(x, t) = K(a,t)W(t), necessitates that STF (2, t)| <M 
for a real number M and for all x,t € [a,b]. Therefore from (27) and (29), 


we obtain 
1 b antip 
aa f Goleta) dt 


b= a)M y ntt aa 


(h= 
a ae (30) 


4(n + 1) n 
where we used (n + 1)n"t? > nn™*! = n”*? to get the last inequality. We 


note that 5 ‘ 
Ge aeG- Gas! 
< 14 : 
nr nr n 


Thus from (30), we obtain 


M (b-a\’ b—a\" 
0<|u(e)-un(o)l< F(*=*) (14 °=*) (31) 
4 n n 
It is clear to see that : 
b-— 
lim ( *) =0, (32) 
noo n 
b— nm 
lim (1+ *) = eb, (33) 
noo n 


Considering the limit of the both sides of (31) as n approaches infinity and 
using (32)—(33), yield 


392 Shahsavaran 


lim un(x) = u(z), for every x € [a, }); 
n> co 


that is, u,(#) is uniformly convergent to u(x). 


As a result of Theorem 3 and from (30), it is also clearly seen that 


—a 


ite) ate OG); where k= 2 (34) 


n 


4 Illustrative examples 


In this section, we apply the method proposed in Section 2 to some test 
examples. All numerical calculations are performed by Maple 13. 


4.1 Fredholm integral equation 


Example 1. Consider 


with exact solution u(x) = 2. 


For this example, by solving the nonlinear system (16) for n = 4, we have 
Wao = 0, War — 0.0625, W492 = 0.25, WVa3 = 0.5625, Waa — il, 


Also the nodes x; and the weights w; of the Newton—Cotes quadrature rule 
for the same n are 


1 i 3 

ti 0, t t t t 1 

0 > 1 4? 2 2? 3 4’ 4 , 

7 16 2 16 is 
Ww 5 W > W > —,W4 = —. 
Og age Be Be 96 


Substituting the values of x;, w; and V4, for i =0,...,4 into (18), we have 
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4 
ua(x) =v (a) + So wits(a, ti) Ua, 
i=0 
1 
aint ee x 0) x 0 
16,1 1 
2 


“is 
16 


* 45 
7 


* 90 


=2, 


1 
=) x 0.25 
x 5) x 


1 

~ 2x 

2 

1 3 

at x D x 0.5625 
1 

pre Si 

2 


x1)x1 


which is the exact solution of Example 1 


Example 2. [9] Consider 
1 
u(x) =er+1 -f (at+tye"Odt, O<2<1, 
0 


with the exact solution u(x) = a. 


The absolute error |u,,(a) — u(x)| of Example 2 for n = 2, 4,6 is shown in 
Table 1. 


Table 1: Absolute error of Example 2 


x n=2 | n=4 | n=6 | method of [9] (N =7) 
0.0 | 1.7e—3 | 3.8e—6 | 6.0e—9 2.5e—6 
0.2 | 1.5e-3 | 3.2e—6 | 5.2e—9 7.3e—6 
0.4 | 1.2e-3 | 2.6e—6 | 4.4e—9 7.9e—6 
0.6 | 1.0e—-3 | 2.le—6 | 3.6e—9 2.5e—6 
0.8 | 7.7e—4 | 1.5e-6 | 2.8e—9 3.9e—6 
1.0 | 5.3e—4 | 9.6e—7 | 2.0e—9 2.6e—6 


Example 3. [1] Consider 


with exact solution u(a) = e*. 


The absolute error |u,(a) — u(x)| of Example 3 for n = 2, 4,6 is shown in 
Table 2. 
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Table 2: Absolute error of Example 3 


x n=2 | n=4 | n=6 | method of [1] (NV = 16) 
0.1 | 1.0e—4 | 1.5e—7 0.0 4.3e—4 
0.2 | 1l.le—4 | 1.7e—7 | 1.0e—9 3.3e—7 
0.3 | 1.2e—4 | 1.8e—7 0.0 7.9e—4 
0.4 | 1.4e—4 | 2.0e—7 0.0 2.9e—4 
0.5 | 1.5e—4 | 2.3e—7 | 1.0e—9 1.2e—3 
0.6 | 1.7e—4 | 2.5e—7 0.0 7.le—4 
0.7 | 1.8e—4 | 2.8e—7 | 1.0e—9 5.4e—7 
0.8 | 2.0e—4 | 3.le—7 | 1.0e—9 1.3e—3 
0.9 | 2.3e—4 | 3.4e—7 | 1.0e—9 4.8e—4 


Example 4. [12] In this example, the proposed method in Section 2 is used to 
solve an integral equation reformulation of the nonlinear two-point boundary 
value problem 


u (t) —exp(u(t))=0, t£€(0,1), u(0)=u(1) =0. 


This problem has the unique solution 


u(t) = —In(2) + 2In <=“) 


= /2, and may be reformulated as the 


Cc 


cos(¢) 


where c is the only solution of 


integral equation 


u(x) = | k(x, t) exp(u(t))dt, a € [0,1], 


where 


eae Hee t<a, 


—a(1—-t), t>za. 


The uniform norm ||un — u|| = sup{|un(t) — u(t), t € [0, 1]} of Example 
4 for n = 5,9 is shown in Table 3. 


Table 3: Numerical results of Example 4 
n | ||U, —ul| (presented method) | ||u, — ul] (method of [12]) 
5.94e—3 5.19e—4 
9 2.19e—3 1.28e—4 


on 
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4.2 Volterra integral equation 


Example 5. [5] Consider 


with the exact solution u(#) = e7*. 


For this example, we define e, (2) = |un(x) — u(x)|. The maximum norm 
llen|loo = max{le,(2;)|, v7; = .1*%,i = 0,1,...,10} for n = 1,3,5,7,9 is 
presented in Table 4. 


Table 4: Numerical results of Example 5 


n | |len|loo (presented method) | |{en|]oo (method of [5]) 
1 1.174e—1 6.282e—2 
3 5.235e—4 1.001le—3 
5 3.829e—6 8.294e—6 
7 2.020e—8 3.913e—8 
9 1.000e—10 1.163e—10 


Example 6. [18] Consider 
u(x) =1+sin? 2 — 3 | sin(x — t)(u(t))*dt, O<a2<1, 
0 


with exact solution u(x) = cosa. 


The absolute error |u,(x) — u(a)| of Example 6 for n = 3,5,7,9 is shown 
in Table 5. 


Table 5: Absolute error of Example 6 
x n= n=5 | n=7 n=9 | method of [18] (m = 60) 
0.2 | 4.6e—3 | 8.4e—5 | 8.4e-7 | 6.3e—9 0.0 
0.4 | 3.8e—3 | 2.2e—5 | 4.2e—-7 | 3.9e—9 1.0e—5 
0.6 | 2.le—3 | 1.0e—6 | 1.le—7 | 5.0e—10 1.0e—5 
0.8 | 3.4e—3 | 6.3e—5 | 6.le—7 | 3.9e—9 2.0e—5 
1.0 | 1.le—3 | 2.2e—5 | 2.9e—7 | 2.le—9 1.0e—5 


Example 7 (Constructed by author). Consider 


u(z) =e *—e*(c+1)+ | et ult)dt, -1<a<l. (35) 
<1 
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To calculate the error in the interval [—1,1], we define the error function 


€n(x) as 
xz 


€n(£) = Un(a) — v(x) — / K(x, t)Up(t)dt. 
-1 
Actually, on the right-hand side of the above equation, we put the ap- 
proximated solution uy (x) instead of the exact solution u() for (4). Now the 
absolute error |e,(a)| of Example 7 for n = 2,3,5,7,9 and some x € [—1, 1] 
is shown in Table 6. 


Table 6: Absolute error of Example 7 

x €2(x) €3 (2x) €5 (2x) €7(x) €9(x) 
-1.0 0.0 0.0 0.0 0.0 0.0 
-0.8 | 2.0e—3 | 3.4e—4 | 2.3e—4 | 5.7e—6 | 1.0e—9 
-0.6 | 1.0e—2 | 1.2e—3 | 4.5e—4 | 6.9e—-6 | 1.0e—9 
-0.4 | 2.8e—2 | 2.2e—3 | 4.4e—4 | 3.6e—-7 | 1.2e—9 
-0.2 | 5.2e—2 | 2.5e—3 | 4.2e—4 | 1.2e—5 | 9.0e—10 
0.0 | 7.4e—-2 | 1.6e—-3 | 6.4e—4 | 1.9e—-5 | 3.0e—9 
0.2 | 7.38e-2 | 5.9e—5 | 9.4e—4 | 3.2e—-5 | 6.le—9 
0.4 | 2.0e—2 | 3.7e—4 | 9.0e—4 | 1.8e—5 | 9.2e—9 
0.6 | 1.0e—1 | 1.5e—3 | 7.le—4 | 6.9e-6 | 2.le—8 
0.8 | 2.9e—1 | 2.2e—4 | 1.4e—3 | 3.8e—5 | 6.9e—8 
1.0 | 4.4e—-1 | 1.3e—-3 | 1.4e-3 | 5.le—5 | 2.2e—7 


Example 8. Consider 


vw Qa 5x3 


u(x) = 1 3 


x? +1 +/ (at +1)(u(t))?dt, O<a<1. 
0 


Similar to Example 7, to calculate the error in the interval [0,1], we define 
the error function e,,(2) as 


€n(L) = Un(x) — v(x) — / K(x, t)(Un (t))7 dt. 
0 
Now the absolute error |e,(2)| of Example 8 for n = 2,3,5,7 and some 
x € [0,1] is shown in Table 7. 
5 Conclusion 


In this work, the Newton—Cotes quadrature rule together with the Lagrange 
interpolation were used to transform Fredholm and Volterra integral equa- 
tions to a system of algebraic equations. As shown in Section 4, via some 
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Table 7: Absolute error of Example 8 


x €2(x) €3(2x) €5(x) €7(x) 
0.0 0.0 0.0 0.0 0.0 
0.1 | 4.7e—4 | 9.5e—3 0.0 0.0 


0.2 | 1l.le—3 | 2.7e—2 | 1.0e—10 | 1.0e—10 
0.3 | 1.4e—3 | 4.0e—2 | 1.0e—10 | 3.0e—10 
0.4 | 1l.le—3 | 4.1le—2 | 1.0e—10 | 2.0e—10 
0.5 | 4.4e—4 | 3.le—2 | 1.0e—10 0.0 

0.6 | 9.4e—4 | 1.9e—2 0.0 2.0e—9 
0.7 | 3.4e—3 | 1.7e—2 0.0 8.0e—9 
0.8 | 7.4e—3 | 3.2e—2 | 1.0e—9 1.0e—8 
0.9 | 1.2e—2 | 5.le—2 | 1.0e-9 | 4.0e—8 
1.0 | 1.6e—2 | 5.7e—2 | 3.0e—9 7.7e—8 


test examples, as n increases, the error decreases. The calculated errors in 
the test examples were compatible with the presented error bound in (30). 
It was stated that a high accuracy is achieved even by using a small number 
of n. Also the method can be extended to solve a system of such equations. 
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